﻿/***********************************************************************
 * CLR 版本: 4.0.30319.225
 * 工程名称: $projectname$
 * 命名空间: Toulr.Util.Geography
 * 文件名称: GeoUtils
 * 描    述: 
 * 作    者: 18901309129@189.cn
 * 登录用户: hyyu
 * 创建时间: 2011-6-14 18:56:42
 * 创建年份: 2011
 * Update History: 
 ***********************************************************************/

using System;

namespace Toulr.Util.Geography
{
    public class GeoUtils
    {
        /// <summary>
        /// 地球极半径，单位为m
        /// </summary>
        public const double EarthPolarRadius = 6356725; //为了跟后台计算统一,据wiki为6356752.3

        /// <summary>
        /// 地球平均半径，单位为m
        /// </summary>
        public const double EarthMeanRadius = 6371009;

        /// <summary>
        /// 给定地面两个点的经纬度，计算其球面距离(大圆上的劣弧)
        /// </summary>
        /// <param name="p1"></param>
        /// <param name="p2"></param>
        /// <returns></returns>
        public static double GetDistanceMeters(GeoPoint p1, GeoPoint p2)
        {
            return GetDistanceMeters(p1.Latitude, p1.Longitude,
                                     p2.Latitude, p2.Longitude);
        }

        /// <summary>
        /// 给定地面两个点的经纬度，计算其球面距离(大圆上的劣弧)，单位是米
        /// </summary>
        /// <param name="lat1">第一个点的纬度</param>
        /// <param name="lng1">第一个点的经度</param>
        /// <param name="lat2">第二个点的纬度</param>
        /// <param name="lng2">第二个点的经度</param>
        /// <returns></returns>
        public static double GetDistanceMeters(double lat1, double lng1,
            double lat2, double lng2)
        {
            double radLat1 = Rad(lat1);
            double radLat2 = Rad(lat2);
            double a = radLat1 - radLat2;
            double b = Rad(lng1) - Rad(lng2);

            double s = 2 * Math.Asin(Math.Sqrt(Math.Pow(Math.Sin(a / 2), 2) +
                Math.Cos(radLat1) * Math.Cos(radLat2) * Math.Pow(Math.Sin(b / 2), 2)));
            s = s * EarthPolarRadius;
            return s;
        }

        /// <summary>
        /// 计算在指定的point附近，length米的距离相当于多少经度，单位为度
        /// </summary>
        /// <param name="point"></param>
        /// <param name="length"></param>
        /// <returns></returns>
        public static double GetLongitudeDelta(GeoPoint point, int length)
        {
            double lngRadius = EarthPolarRadius * Math.Cos(Rad(point.Latitude));
            return 180 * length / (Math.PI * lngRadius);
        }

        /// <summary>
        /// 计算在指定的point附近，length米的距离相当于多少纬度，单位为度
        /// </summary>
        /// <param name="point"></param>
        /// <param name="length"></param>
        /// <returns></returns>
        public static double GetLatitudeDelta(GeoPoint point, int length)
        {
            return 180 * length / (Math.PI * EarthPolarRadius);
        }

        /// <summary>
        /// 将角度从度单位转换为弧度单位
        /// </summary>
        /// <param name="d"></param>
        /// <returns></returns>
        public static double Rad(double d)
        {
            return d * Math.PI / 180.0;
        }

        /// <summary>
        /// 判断一点是否在多边形内部
        /// </summary>
        /// <param name="pt">目标点</param>
        /// <param name="polygonPoints">多边形端点数组</param>
        /// <returns></returns>
        public static bool IsPointInsidePolygon(GeoPoint pt,
                                                  GeoPoint[] polygonPoints)
        {
            // 多边形最少3个边
            if (polygonPoints == null || polygonPoints.Length < 3)
            {
                return false;
            }

            // 横向交点个数
            double dCorssCountX = 0;

            // 设定原点，为第一个端点的坐标
            GeoPoint gCenter = polygonPoints[0];

            // 目标点经纬度转换为直角坐标
            PointDouble dTargetPoint = GetTriangleXAndY(gCenter, pt);

            // 直角坐标的端点
            PointDouble[] dSidePoint = new PointDouble[polygonPoints.Length];

            // 计算每个端点的直角坐标系坐标
            for (int nIndex = 0; nIndex < polygonPoints.Length; nIndex++)
            {
                dSidePoint[nIndex] = GetTriangleXAndY(gCenter, polygonPoints[nIndex]);
            }

            // 判断目标点是否在边上，若在边上，直接返回true
            for (int nIndex = 0; nIndex < dSidePoint.Length; nIndex++)
            {
                // 最后一次循环
                if (nIndex == dSidePoint.Length - 1)
                {
                    if (IsPointInLine(dSidePoint[nIndex], dSidePoint[0], dTargetPoint))
                    {
                        return true;
                    }
                }
                else
                {
                    if (IsPointInLine(dSidePoint[nIndex], dSidePoint[nIndex + 1], dTargetPoint))
                    {
                        return true;
                    }
                }
            }

            // 分别计算目标点和多边形边的交点个数
            for (int nIndex = 0; nIndex < dSidePoint.Length; nIndex++)
            {
                // 横向每次计算的交点数
                double dCrossCountTempX;

                // 最后一次循环
                if (nIndex == dSidePoint.Length - 1)
                {
                    dCrossCountTempX = FindCrossPointCountOfOneSide(
                        dSidePoint[nIndex], dSidePoint[0], dTargetPoint, nIndex, dSidePoint);
                }
                else
                {
                    dCrossCountTempX = FindCrossPointCountOfOneSide(
                        dSidePoint[nIndex], dSidePoint[nIndex + 1], dTargetPoint, nIndex, dSidePoint);
                }
                dCorssCountX += dCrossCountTempX;
            }

            // 转换为int值
            int nCrossCountX = Convert.ToInt32(dCorssCountX);

            return nCrossCountX % 2 != 0;
        }

        /// <summary>
        /// 根据地理经纬度换算直角坐标
        /// </summary>
        /// <param name="center">经纬度原点</param>
        /// <param name="target">要换算的点</param>
        /// <returns>目标点的直角坐标</returns>
        private static PointDouble GetTriangleXAndY(GeoPoint center,
                                                    GeoPoint target)
        {
            // 参数错误，返回null
            if (center.IsEmpty)
            {
                return PointDouble.Empty;
            }

            double x = 2 * Math.PI * EarthPolarRadius *
                       Math.Cos(center.Latitude * Math.PI / 180) *
                       (target.Longitude - center.Longitude) / 360;
            double y = 2 * Math.PI * EarthPolarRadius *
                       (target.Latitude - center.Latitude) / 360;
            PointDouble dPoint = new PointDouble(x, y);
            return dPoint;
        }

        /// <summary>
        /// 计算点和边的交点数（自点向左作射线）
        /// </summary>
        /// <param name="sidestart">边的起点</param>
        /// <param name="sideend">边的终点</param>
        /// <param name="point">要判断的点</param>
        /// <param name="index">索引</param>
        /// <param name="pointlist">端点集合</param>
        /// <returns>交点数量</returns>
        private static float FindCrossPointCountOfOneSide(
            PointDouble sidestart, PointDouble sideend, PointDouble point,
            int index, PointDouble[] pointlist)
        {

            // 参数错误，返回false
            if (sidestart.IsEmpty || sideend.IsEmpty || point.IsEmpty)
            {
                return 0.0f;
            }

            // 边竖直
            if (sidestart.X == sideend.X)
            {
                double dVerticalMaxY;
                double dVerticalMinY;
                if (sidestart.Y > sideend.Y)
                {
                    dVerticalMaxY = sidestart.Y;
                    dVerticalMinY = sideend.Y;
                }
                else
                {
                    dVerticalMaxY = sideend.Y;
                    dVerticalMinY = sidestart.Y;
                }
                // 点的横坐标大于边的横坐标
                if (point.X > sidestart.X)
                {
                    // 点的纵坐标在两端点纵坐标范围内，一个交点
                    if (point.Y > dVerticalMinY && point.Y < dVerticalMaxY)
                    {
                        return 1.0f;
                    }
                    // 点在边的起点或终点上
                    else if (point.Y == dVerticalMinY || point.Y == dVerticalMaxY)
                    {
                        return CalcCrossCount(sidestart, sideend, point, index, pointlist);
                    }
                    // 点不在边纵坐标范围内，无交点
                    else
                    {
                        return 0.0f;
                    }
                }
                // 点的横坐标小于边的横坐标
                else
                {
                    return 0.0f;
                }
            }

            // 边水平
            if (sidestart.Y == sideend.Y)
            {
                return 0.0f;
            }

            // 三个点的坐标
            double x1 = sidestart.X;
            double y1 = sidestart.Y;
            double x2 = sideend.X;
            double y2 = sideend.Y;
            double x3 = point.X;
            double y3 = point.Y;
            // 求交点
            double crossPointY = y3;
            double crossPointX = ((x2 - x1) * (y3 - y1)) / (y2 - y1) + x1;
            PointDouble dCrossPoint = new PointDouble(crossPointX, crossPointY);

            // 点在交点左边，无交点
            if (dCrossPoint.X > point.X)
            {
                return 0.0f;
            }
            else
            {
                // 交点在线段上，1个交点
                if (CheckPointInLine(sidestart, sideend, dCrossPoint))
                {
                    // 交点是端点时
                    if ((dCrossPoint.X == sidestart.X && dCrossPoint.Y == sidestart.Y) ||
                        (dCrossPoint.X == sideend.X && dCrossPoint.Y == sideend.Y))
                    {
                        return CalcCrossCount(sidestart, sideend, point, index, pointlist);
                    }
                    // 交点不是端点，计1.0
                    else
                    {
                        return 1.0f;
                    }
                }
                // 交点不在线段上，无交点
                else
                {
                    return 0.0f;
                }
            }
        }

        /// <summary>
        /// 两边交点也是射线交点时，数量计算
        /// </summary>
        /// <param name="sidestart">边起点</param>
        /// <param name="sideend">边终点</param>
        /// <param name="point">目标点</param>
        /// <param name="index">边起点索引</param>
        /// <param name="pointlist">端点集合</param>
        /// <returns>交点计数</returns>
        private static float CalcCrossCount(PointDouble sidestart, PointDouble sideend,
                                             PointDouble point, int index, PointDouble[] pointlist)
        {
            // 交于起点
            if (point.Y == sidestart.Y)
            {
                // 终点比交点大
                if (sideend.Y > point.Y)
                {
                    // 不是第一条
                    if (index > 0)
                    {
                        // 前边起点比交点大
                        if (pointlist[index - 1].Y >= point.Y)
                        {
                            return 1.0f;
                        }
                        else
                        {
                            return 0.5f;
                        }
                    }
                    // 第一条的时候，和最后一条计算
                    else
                    {
                        // 前边起点比交点大
                        if (pointlist[pointlist.Length - 1].Y >= point.Y)
                        {
                            return 1.0f;
                        }
                        else
                        {
                            return 0.5f;
                        }
                    }
                }
                else
                {
                    // 不是第一条
                    if (index > 0)
                    {
                        // 前边起点比交点小
                        if (pointlist[index - 1].Y <= point.Y)
                        {
                            return 1.0f;
                        }
                        else
                        {
                            return 0.5f;
                        }
                    }
                    // 第一条的时候，和最后一条计算
                    else
                    {
                        // 前边起点比交点大
                        if (pointlist[pointlist.Length - 1].Y <= point.Y)
                        {
                            return 1.0f;
                        }
                        else
                        {
                            return 0.5f;
                        }
                    }
                }
            }
            // 交于终点
            else
            {
                // 起点比交点大
                if (sidestart.Y > point.Y)
                {
                    // 不是最后2条边
                    if (index < pointlist.Length - 2)
                    {
                        // 后边起点比交点大
                        if (pointlist[index + 2].Y >= point.Y)
                        {
                            return 1.0f;
                        }
                        else
                        {
                            return 0.5f;
                        }
                    }
                    // 倒数第二条边，和第一条计算
                    else if (index == pointlist.Length - 2)
                    {
                        // 后边起点比交点大
                        if (pointlist[0].Y >= point.Y)
                        {
                            return 1.0f;
                        }
                        else
                        {
                            return 0.5f;
                        }
                    }
                    // 最后一条边，和第二条计算
                    else
                    {
                        // 后边起点比交点大
                        if (pointlist[1].Y >= point.Y)
                        {
                            return 1.0f;
                        }
                        else
                        {
                            return 0.5f;
                        }
                    }
                }
                else
                {
                    // 不是最后2条边
                    if (index < pointlist.Length - 2)
                    {
                        // 后边起点比交点小
                        if (pointlist[index + 2].Y <= point.Y)
                        {
                            return 1.0f;
                        }
                        else
                        {
                            return 0.5f;
                        }
                    }
                    // 倒数第二的时候，和第一条计算
                    else if (index == pointlist.Length - 2)
                    {
                        // 后边起点比交点大
                        if (pointlist[0].Y <= point.Y)
                        {
                            return 1.0f;
                        }
                        else
                        {
                            return 0.5f;
                        }
                    }
                    // 最后一条边，和第二条计算
                    else
                    {
                        // 后边起点比交点大
                        if (pointlist[1].Y <= point.Y)
                        {
                            return 1.0f;
                        }
                        else
                        {
                            return 0.5f;
                        }
                    }
                }
            }
        }

        /// <summary>
        /// 判断点是否在线段内
        /// </summary>
        /// <param name="startline">线段起点</param>
        /// <param name="endline">线段终点</param>
        /// <param name="crosspoint">要判断的点</param>
        /// <returns>是否在线段上</returns>
        private static bool CheckPointInLine(PointDouble startline, PointDouble endline,
            PointDouble crosspoint)
        {
            // 参数错误，返回false
            if (startline.IsEmpty || endline.IsEmpty || crosspoint.IsEmpty)
            {
                return false;
            }

            // 点的大小
            double dMinX;
            double dMaxX;
            double dMinY;
            double dMaxY;

            // X点的比较
            if (startline.X > endline.X)
            {
                dMinX = endline.X;
                dMaxX = startline.X;
            }
            else
            {
                dMinX = startline.X;
                dMaxX = endline.X;
            }
            // Y点的比较
            if (startline.Y > endline.Y)
            {
                dMinY = endline.Y;
                dMaxY = startline.Y;
            }
            else
            {
                dMinY = startline.Y;
                dMaxY = endline.Y;
            }
            // 点在不在线段上的判断
            if (crosspoint.X > dMaxX || crosspoint.X < dMinX || crosspoint.Y > dMaxY || crosspoint.Y < dMinY)
            {
                return false;
            }
            else
            {
                return true;
            }
        }

        /// <summary>
        /// 判断P点是否在P1P2这条线段上
        /// </summary>
        /// <param name="p1">线段起点</param>
        /// <param name="p2">线段终点</param>
        /// <param name="p">目标点</param>
        /// <returns>点是否在线段上的结果</returns>
        private static bool IsPointInLine(PointDouble p1, PointDouble p2, PointDouble p)
        {
            // P点是边的两端点之一，返回true
            if (p.X == p1.X && p.Y == p1.Y || p.X == p2.X && p.Y == p2.Y)
            {
                return true;
            }

            double x1 = p1.X;
            double y1 = p1.Y;
            double x2 = p2.X;
            double y2 = p2.Y;
            double x3 = p.X;
            double y3 = p.Y;

            double dXMax;
            double dXMin;
            double dYMax;
            double dYMin;

            if (x1 > x2)
            {
                dXMax = x1;
                dXMin = x2;
            }
            else
            {
                dXMin = x1;
                dXMax = x2;
            }
            if (y1 > y2)
            {
                dYMax = y1;
                dYMin = y2;
            }
            else
            {
                dYMin = y1;
                dYMax = y2;
            }

            // 线段竖直
            if (x1 == x2)
            {
                // 在直线上
                if (x3 == x1)
                {
                    // 在线段上
                    if (y3 >= dYMin && y3 <= dYMax)
                    {
                        return true;
                    }
                    else
                    {
                        return false;
                    }
                }
                else
                {
                    return false;
                }
            }

            // 线段水平
            if (y1 == y2)
            {
                // 在直线上
                if (y3 == y1)
                {
                    // 在线段上
                    if (x3 >= dXMin && x3 <= dXMax)
                    {
                        return true;
                    }
                    else
                    {
                        return false;
                    }
                }
                else
                {
                    return false;
                }
            }

            // 斜率
            double dKPP1 = (y3 - y1) / (x3 - x1);
            double dKPP2 = (y3 - y2) / (x3 - x2);

            // 斜率相同，则在同一直线上
            if (dKPP1 == dKPP2)
            {
                // 在线段上
                if (x3 >= dXMin && x3 <= dXMax && y3 >= dYMin && y3 <= dYMax)
                {
                    return true;
                }
                else
                {
                    return false;
                }
            }
            else
            {
                // 斜率不同，不在一条直线上
                return false;
            }
        }
    }
}
